#include <iostream>
#include <Eigen/Dense> // 包含完整的Eigen库
#define MAX_ITER 400

using namespace std;

int main()
{
    // 创建一个3x3的矩阵
    Eigen::Matrix3d A;
    A << 1, -1, 0, -1, 2, -1, 0, -1, 1;
    Eigen::Vector3d b;
    b << -1, -1, 2;

    const double atol = 1e-8;
    const double step_size = 0.1;  // 迭代步长
    const float eps[6] = {1, 1e-1, 1e-2, 1e-3, 1e-4, 0};
    Eigen::Vector3d x = Eigen::Vector3d::Zero();  // 初始化x0
    Eigen::Vector3d g = Eigen::Vector3d::Zero();


    // 对角矩阵
    Eigen::Matrix3d I = Eigen::Matrix3d::Identity();
    I.diagonal() << 1, 1, 1; // 设置对角线元素

    // 循环不同的epsilon 比较迭代步数
    for (int i=0; i<6; i++){
        printf("\neps = %e\n", eps[i]);
        auto A_eps = A + eps[i] * I;
        // 梯度下降算法
        int iter = 0;  // 迭代次数
        // norm默认计算l2范数
        Eigen::Vector3d r = Eigen::Vector3d::Ones(); // 初始化残差r
        while (iter < MAX_ITER and r.norm() > atol)
        {
            r = A_eps * x - b;
            g = A_eps.transpose()*r;
            x -= step_size * g;
            iter++;
        }
        printf("Iter=%i, r=%e\n", iter, r.norm());
        cout << "x =\n" <<  x << endl;
    }
    return 0;
}